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Abstract We study the transition from regular to chaotic motion in a prolate elliptical 
galaxy dynamical model with a bulge and a dense nucleus. Our numerical investigation 
shows that stars with angular momentum less than or equal to a critical value L^c, 
moving near the galactic plane, are scattered to higher z, when reaching the central region 
of the galaxy, thus displaying chaotic motion. An inverse square law relationship was 
found to exist between the radius of the bulge and the critical value L^c of the angular 
momentum. On the other hand, a linear relationship exists between the mass of the nucleus 
and Lzc- The numerically obtained results are explained using theoretical arguments. Our 
study shows that there are connections between regular or chaotic motion and the physical 
parameters of the system, such as the star's angular momentum and mass, the scale length 
of the nucleus and the radius of the bulge. The results are compared with the outcomes of 
previous work. 
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1 INTRODUCTION 

Today astronomers believe that the intrinsic shapes of most elliptical galaxies are either oblate or prolate. 



It is also true that spherical galaxies are very rare but triaxial elliptical galaxies do exist (see Statler 1994 
|Ryden|1996l|Alani & Ry den 2002 , [ Vincent & Ryden|2005| l. On the other hand, there are observational 
data indicating the presence of a black hole or a dense massive nucleus in the central parts of elliptical 
galaxies (Statler et al.,2004j . 

All of the above gives us the opportunity to construct a dynamical model for an elliptical galaxy 
hosting a dense nucleus, in order to use it for the study of the global properties of motion in these stellar 
systems. In order to describe the motion in a prolate elliptical galaxy we use the potential 

n^fln{r' + az' + cl)-—M^=V, + V^. (1) 

Our model consists of two components. The first component describes a prolate elUptical galaxy while 
the second is the potential of a dense spherical nucleus. Here (r, z) are the usual cylindrical coordinates, 
vq is used for consistency of galactic units, 0.2 < a < 1 is the flattening parameter, and Cb is the radius 
of the bulge component. Furthermore, M„ is the mass and is the scale length of the nucleus. 

In an earlier paper ( Caranicolas & Innanen|199T hereafter Paper I) we have studied the transition 



from regular to chaotic motion in a disk galaxy model with a dense nucleus. There we found that stars, 
moving in the (r, z) plane with values of angular momentum less than or equal to a critical value 
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Lzc, are scattered to higher z upon encountering the dense nucleus, thus displaying chaotic motion. We 
also found relationships connecting the physical parameters of the system with chaos. In the present 
paper, we shall focus our study on the transition from regular to chaotic motion in model (1), which 
describes a prolate elliptical galaxy hosting a dense nucleus. Our aim is: (i) to look for relationships 
between the physical parameters and chaos, (ii) to explain the numerically found relationships using a 
semi-theoretical approach and (iii) to compare the present results with the outcomes found in Paper I. 

In this research we use the well known system of galactic units, where the unit of length is 1 kpc, the 
unit of mass is 2.325 x 10^ Mq and the unit of time is 0.97748 x 10^ yr. The velocity and the angular 
velocity units are 10 km s^^ and 100 km s^^ kpc^^ respectively, while G is equal to unity. The energy 
unit (per unit mass) is 100 km^ s^^. In the above units we use the values: — 10 while < < 150, 
0.75 < Cb < 1.25 and 0.1 < < 0.25. Since the total potential Vt — Vt (r, z) is axially symmetric 
and the component of the angular momentum is conserved, we use the effective potential 

V,s^ ^+Vtir,z), (2) 
2r"^ 

in order to study the motion in the meridian (r — z) plane. The equations of motion are 



Pr = ^ , Pz= 5 , (3) 

or oz 
and the corresponding Hamiltonian is written as 

^=^(P'+P')+Kff(r,z)=i?, (4) 

where p^ and pz are the momenta per unit mass conjugate to r and z respectively, while E is the 
numerical value of the Hamiltonian. Equation (4) is an integral of motion, which indicates that the 
total energy of the particle is conserved. Orbit calculations are based on the numerical integration of 
the equations of motion (3), which were made using a Bulirsh-Stoer routine, to double precision. The 
accuracy of the calculations was checked by the constancy of the energy integral, which was conserved 
up to the twelfth significant figure. 

The paper is organized as follows. In Section 2 we present numerical results for the potential of 
a non-active galaxy, that is when = 0. Furthermore, the numerical relationship is explained using 
some semi-theoretical arguments. In Section 3, we investigate numerically the case when the dense 
nucleus is present. Some semi-theoretical arguments are also presented, in order to explain the numerical 
outcomes. In Section 4 a comparison with earlier work is given and discussion and conclusions of this 
research are presented. 



2 RESULTS WHEN THE MASSIVE DENSE NUCLEUS IS NOT PRESENT 

In this Section we shall study the behavior of orbits when the dense massive nucleus is not present, that 
is when Mn ~ 0. Figure[T]shows the numerical relationship between the critical value Lzc of the angular 
momentum and the radius of the bulge Cb- Orbits were started near rg — r^y^^^, with zq = Pro — 0, while 
the value of pzo is always found from the energy integral (4). The value of rmax is the maximal root of 
equation 

^^+\\n{r' + cl)=E, (5) 

which was found numerically. Dots represent the numerical values while the solid line is the best fit 
which is an inverse square law, represented by the equation 
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Fig. 1 Plot of vs. Cb. The values of the parameters are given in the 
text. 
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Fig. 2 (r.pr) phase plane when the dense nucleus is absent, (a) Cb = 0.75 and (b) Cb = 1.25. 
The values of all other parameters are given in the text. 



The value of a is 0.2, while the value of i? is 231. Orbits with values of the parameters on and below 
the line are chaotic, while orbits on the upper side of the line are regular 

Figure |2] shows the Poincare {r,pr), z — 0,pz > phase plane when a = 0.2, — I, E — 213, 
while the value of Cf, is 0.75 in Figure |2| a) and 1.25 in Figure |2|b). As one can see, the pattern is 
similar in both figures with regions of regular and chaotic motion. A more detailed inspection shows 
that the extent of the chaotic zone is larger in the case when Cf, = 0.75, which is when we have a 
galaxy with a denser bulge. Figure |3] shows two orbits when a — 0.2, — 0.9, E — 231. The value 
of the angular momentum in the orbit shown in Figure [3j a) is ~ 40 and the initial conditions are: 
To = 9.1, zq = Pro — 0, while the value of pzo is always found from the energy integral (4). Note 
that the orbit is regular and stays very close to the galactic plane. On the contrary, the orbit shown in 
Figure[3]^b) has a value of — 8 and initial conditions: tq = 10.0, zq — p^o — 0. The orbit is chaotic 
and it is scattered off the galactic plane, displaying high values of z. Both orbits were calculated for a 
time period of 200 time units. 

Let us now begin to use semi-analytical arguments in order to explain the numerically found rela- 
tionship of Figure [T| The lines of arguments are similar to those used in Paper I. As the test particle 
approaches very close to the center of a galaxy, there is a change in its momentum in the z direction 
given by the equation 

mAp, - {F,)At. (7) 
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Fig, 3 Orbits when the dense nucleus is not present, (a) A regular orbit which stays near the 
galactic plane, (b) A chaotic orbit which is scattered into high values of z. See tiie text for 
details. 



Here m is the mass of the test particle, (Fj) is the total average force acting in the z direction and 
At is the duration of the encounter. Empirical evidence shows that the test particle's rise proceeds 
cumulatively in each case, and increases a little more with each successive pass from the center rather 
than with a single "violent" encounter. It is observed that the test particle gains considerable height after 
n(n > 1) crossings, when the total change in the momentum in the z direction is on the order of mu^, 
where is the tangential velocity of the test particle near the center at a distance r = (ro) — (^o) 1- 
Therefore we write 



m 



^Ap,,«(F,)^Ati. (8) 



Setting 



m 



in Equation (8) we find 



mL, 



(ro) ' 

To, (9) 

, , {F,)T,. (10) 
The force acting in the z direction for a test particle of unit mass {m = 1) is 

r-^ + az^ + 

Remember that r = (ro) — (^o) *C 1, therefore (rg) ~ {zq) <C c^. Keeping only the linear terms in (r) 
and (z) in Equation (11) and taking the absolute value of the force, we find from relationship (10) 
that 

i;ga((ro))^^ 



''ZC 



(12) 



Since Tc was observed to be the same when 0.75 < Cb < 1.25, we can set k\ = Wgo; ((^o))^ and 
obtain ^ 

^zc ~ (13) 
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Fig. 4 Plot of Mn vs. L^c- The values of the parameters are given in 
the text. 



20 




Fig. 5 (r, pr) phase plane when the dense nucleus is present, (a) Mn = 20 and (b) = 100. 
The values of all other parameters are given in the text. 

which explains the numerical relationship of Figure [T] The authors would like to make it clear that 
relation (13) does not reproduce the numerical results shown in Figure[T] but only shows the form of the 
relationship. 

3 RESULTS WHEN THE MASSIVE NUCLEUS IS PRESENT 

We investigate the case when we have an active galaxy, that is when the dense massive nucleus is present. 
Figure |4] shows a numerical relationship between the critical value L^c of the angular momentum and 
the mass of the nucleus for two values of Cn. The procedure to obtain the results shown in Figure |4] 
is similar to that followed in Figure[T] The value of a is 0.5, Cb — 1.2, while the value of E is 227. Here 
we see a straight line. Orbits with values of the parameters on the left side of the line, including the line, 
are chaotic, while orbits on the right side of the line are regular. It is interesting to note that the extent 
of the chaotic region is larger when the value of Cn is smaller, which is when we have a denser nucleus. 

Figure |5ja)-(b) shows the Poincare {r,pr),z = 0,Pz > phase plane when: a — 0.2, Lz = 
10, E = 230, Cb = 1.2; the value of is 20 in Figure [Sj a) and 100 in Figure [Sjb). In both cases, 
we observe regular regions together with large chaotic regions. Some small islands are also present 
indicating secondary resonances. It is evident that the area covered by chaotic orbits is larger in the 
case of A/n — 100, which is when we have a more massive nucleus. In Figure |6|a)-(b) we can see 
two orbits when a = 0.5, Cb = 1.2, E ~ 227. The value of the angular momentum in the orbit shown 
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Fig. 6 Orbits when the dense nucleus is present, (a) A regular orbit which stays near the 
galactic plane, (b) A chaotic orbit which is scattered into high values of z. See the text for 
details. 
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Fig. 7 Plot of I and jF^^ | vs. z. The values of all the other parameters are given in the text. 



in Figure [6]^a) is = 60, while the value of M„ is 25 and c„ = 0.25. The initial conditions are 
ro — 6.0, zq — pro = 0, while the value of p^o is always found from the energy integral (4). Note 
that the orbit is regular, and stays very close to the galactic plane. On the other hand, the orbit shown in 
Figure|6|b) has a value of Lz = 45, while the value of is 120 and Cn — 0.25. Initial conditions are: 
7'o = 9.768, Zq = pro = 0. The orbit is chaotic and it is scattered off the galactic plane, displaying high 
values of z. Both orbits were calculated for a time period of 200 time units. 

The linear relationship of Figure]?] can be obtained using semi-analytical arguments. As the test 
particle approaches the nucleus it experiences a strong vertical force, due to the presence of the dense 
nucleus. 

Figure ]7] shows a plot of the force as well as the nuclear F^n force as a function of z, near the 
nucleus when ro — 0.1. The values of the parameters are a = 0.5, vq — 10, Cb = 0.8, = 100, = 
0.25. We see that the nuclear force is about 25 times as strong as F^. Because of this strong force, there 
is a change in its momentum in the z direction given by the equation 

mApz = (F,„)At, (14) 

where m is the mass of the test particle, {F-^n) is the average nuclear force acting in the z direction, 
while At is the duration of the encounter Here again, our numerical results show that the test particle 
goes to higher z after n{n > 1) crossings, when the total change in the momentum in the z direction is 
on the order of mu^, where is the tangential velocity of the test particle near the center, at a distance 
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ro — ~ Cn- Therefore we write 

n n 



If we set 



m 

4=1 



(ro) ' 



i=l 

m = 1, 

(3c2) 

A/n « hLzcCn, (17) 

where fc2 = a/27 /Tg. Note that the values of were observed to be about the same when < M„ < 
150. Equation (17) explains the numerically found relationship given in Figure|4] 

4 DISCUSSION 



in Equation (15) we obtain 



Today it is well-known that there are luminous elliptical galaxies hosting dense massive nuclei (see Urry 



et al.|2000| . During the past few years, a large amount of observational data has given a better and more 



detailed picture of these active galaxies ( Barth et al. 2002 Costamante & Ghisellini 2002 ; Falomo et al. 



2002 ; Vagnetti et al. l003] [Falcone et al.||2004[ |Heidt"etaL||2004[ |Bramel et al.|2005[ [Nieppola et al. 
2006,. Zheng et al..2007j . 



In this article we have studied the transition from regular to chaotic motion in a prolate elliptical 
galaxy model. Two cases were studied, which are the case where the nucleus was absent and the case 
where we have an active galaxy hosting a dense massive nucleus. In the first case we found that an 
inverse square law relationship exists between the radius of the bulge of the galaxy and the critical 
angular momentum when all other parameters are kept fixed. This relationship was also reproduced, 
using some semi-analytical arguments. 

On the other hand, it was observed that for larger values of the radius of the bulge, for Cb > 1.5, the 
chaotic regions when the dense nucleus is absent, if any, are negligible when 0.2 < a < 1. This result 
is in agreement with the result found in Paper I, where no chaos was observed when Cb > 1.3 and this 
result was independent of the mass of the bulge. Thus, our results suggest that chaos is observed in disk 
and prolate elliptical galaxies when a dense bulge is present. 

Interesting results are obtained when a dense nucleus is present. In this case the numerically found 
relationship connecting L^c and is linear when all other parameters are kept constant. Furthermore, 
the linear relationship depends on the value of Cn, in such a way that the extent of the chaotic region 
is larger when Cn is smaller, which is when the nucleus has a higher density. The linear relationship 
obtained by the numerical integration of the equations of motion was also found using some semi- 
theoretical arguments, together with numerical evidence. 

The present investigation shows that the results obtained are very similar to those obtained for 
disk galaxies in Paper I. Therefore, we can say that low angular momentum stars approaching a dense 
massive nucleus are deflected to higher z, thus displaying chaotic motion. The similarity of the results 
obtained for different galactic models suggests that it is the strong vertical force near the dense nucleus 
that is responsible for this scattering combined with the star's low angular momentum, which allows the 
star to approach the dense nucleus. 
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Some of the latest discoveries obtained from observational astronomy show that supermassive black 
holes, up to a billion solar masses, inhabit the centers of all massive spheroidal galaxies, independent of 
their visible activity; the supermassive black holes are often quiescent with regard to their own radiation 
but always show dynamical behavior. We believe that, with new data from active galaxies, astronomers 
will be able to construct better dynamical models in the near future, in order to study the properties of 
motion in galaxies and to find interesting relationships connecting chaos with the physical parameters 
of these stellar systems. 
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